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On 1 ABSTRACT 

In the framework of the time-dependent mean field theory, we solve the Vlasov- 
Nordheim equation and study the sensitive dependence of this highly non-linear equa- 
tion on the initial conditions. We find that, when nuclear matter is inside the spinodal 
' region, initially small differences between similar trajectories grow exponentially and 

produce strongly divergent paths resulting in a chaotic evolution. We focus on the fact 
that chaos spontaneously arises without putting any additional noise. We calculate the 
largest Lyapunov exponent for various initial densities. This analysis is performed for 
nuclear matter in two dimensions. It is argued that the mean field chaoticity can cause 
nuclear multifragmentation observed in intermediate energy heavy ion reactions. 
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1 Introduction 

Intermediate energy heavy ion collisions offer a unique possibility to study 
nuclear matter under extreme conditions of density and temperature. In the 
energy range between 20 and 100 MeV/nucleon, the formation of hot nuclei and 
their successive decay, through emission of several intermediate mass fragments 
(IMF : Z > 2) with large fluctuations from event to event, has stimulated a 
huge debate on the main mechanism responsible of this phenomenon, commonly 
named "multifragmentation". 

Both statistical models and dynamical simulations ffl, |2| have been widely used 
providing suggestions on the possible paths followed by the nuclear system along 
the multifragmentation process. Surface instabilities or ring formation || have 
been recently considered as possible steps towards fragment formation. Quite 
surprisingly also the percolation model has been successfully used || , suggesting 
a universal behavior typical of phase transitions and common to different fields. 
The onset of instability in the spinodal zone of the nuclear matter Equation of 
State (EOS) has been studied in ref. |3J], by analyzing the corresponding linear 

1 Talk given by A. Rapisarda. 
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response for imaginary frequency as a function of the wave number and of the 
density, even including a noise of Langevin-type to the collision integral in the 
BUU equation || . However already in the 80s, in one of the first attempts to de- 
scribe multifragmentation Knoll et al. || performed numerical calculations in the 
framework of two-dimensional Time Dependent Hartree Fock (TDHF) equation 
which showed many of the features described by the recent BUU calculations 0. 
In particular by evolving an ensemble of Slater determinants, they found a depen- 
dence on the initial conditions and a spontaneous symmetry breaking which lead 
to large fluctuation from event to event and in the charge or mass distributions. 
This scenario suggests a strong non-linear evolution, as recently suggested also 
in ref.01. 

In order to clarify the mechanism of multifragmentation, we study in this con- 
tribution the dynamics of the Vlasov-Nordheim (VN) equation in the spinodal 
zone. We find unambiguously a strong non-linear evolution which implies a sensi- 
tive dependence on the initial conditions and a spontaneous symmetry breaking 
||. This means that the dynamics in the spinodal region is chaotic, although 
the equation is deterministic: similar trajectories can have a completely differ- 
ent evolution once they enter the spinodal region. We estimate also a Lyapunov 
exponent, which gives the mean rate of divergence between the trajectories in 
density space, in line with the usual analysis of dynamical systems |§. We find 
that two-body collisions do not affect the onset of chaoticity in nuclear many- 
body systems. Finally we draw some conclusions. Further details can be found 
in ref. ||. 



2 Theoretical framework 

2.1 The Vlasov-Nordheim equation 

The numerical analysis is performed by solving the Vlasov-Nordheim equation 

%-{H[f}J} = I[f] (1) 

on a two-dimensional lattice, using the same code of ref.||. In eq. (1) H[f] is the 
self-consistent effective one-body Hamiltonian, while I[f] is the two-body Pauli 
blocked collision integral. The single particle phase space is divided into cells 
and nuclear matter is confined in a large torus with periodic boundary condi- 
tions. The side lengths of the torus are equal to L x = 51 fm and L y = 15 fm. 
In order to achieve sufficient accuracy, very small cells having Sx — ^ fm and 
5p = 40 MeV/c are employed. The time integration step was chosen equal to 0.5 
fm/c. We adopt a simplified Skyrme interaction for the effective one-body field 
which is averaged over the y-direction || and convoluted over a Gaussian func- 
tion whose width is \x = 0.43 fm. This is equivalent to introduce a finite effective 
range. We employ a Fermi momentum at saturation of Pp = 260 MeV/c, which 
gives a density in two dimensions equal to po = 0.55 fm~ 2 . 
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2.2 Definition of the largest Lyapunov exponent 



In the theory of dynamical systems, a measure of chaoticity is expressed in 
terms of Lyapunov exponents || . The largest Lyapunov exponent Aoo can be ex- 
tracted by considering at a given time t the distance d(t) between two trajectories, 
along an unstable direction in phase space, and define 

X(t) = lQgW/c?o) , (2) 

t 

where do = d(0). The average of this quantity along the trajectory gives the 
mean rate of exponential separation of neighbouring trajectories 

Aqo = h m h m A(t) (3) 

t— >oo do—»0 

In numerical applications one has to select a series of small values of do, and 
check that the corresponding values of A(t) stabilize around a definite value of A. 

In order to characterize the system as a whole we have chosen as distance 
between two trajectories the difference in norm between their density profiles 

d(t) = £ \p?\t)-p?\t)\/N c , (4) 

i 

where the index % runs over the N c cells in ordinary space, and p\ , pf 1 are 
the densities in the cell i for the trajectories 1 and 2 respectively The definition 
of eq. (4) is sufficient for the present analysis, although possible differences in 
momentum space are averaged out. It should include the contribution of all the 
unstable modes which dynamically grow up during the evolution. 



3 Results and discussion 

First we analyze the mean field dynamics neglecting two-body collisions, i.e. 
solving the left-hand side of eq.(l). We have initialized nuclear matter with a 
sinusoidal profile along the x-direction, characterized by a wave number k = 
27r(nk/L x ) with = 5 and a small amplitude equal to 1% of the local density. 
The local momentum distribution is assumed to be the one of a Fermi gas at a 
temperature T = 3 MeV. 

In fig. la we show the evolution of the density profile for different times at 
normal density. The initial oscillation along the x coordinate is rapidly damped 
and, as fig. lb shows, the evolution of the Fourier transform of the density profile 
presents only one peak which disappears during the same time scale. Changing 
the initial density by a small amount, for example do = ^ = 10~ 2 , one obtains 
exactly the same behavior and the evolution of the two events is indistinguishable. 
This is typical of the "regular" region of a dynamical system, and indicates the 
stability of the dynamics with respect to small perturbations. 
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On the contrary, if the density p is chosen well inside the spinodal region, 
p = 0.5po in fig-2, the initial symmetry of the density profile, after a short time, 
is completely broken and fluctuations are rapidly enhanced. 

The Fourier transform on the other hand shows the appearance of several 
new modes which where not excited initially. Moreover a small change in the 
initial density, see fig.2bis where p = 0.51p , produces an evolution which is 
completely different both in shape and in amplitude from that of fig.2. This 
is the typical behavior of a dynamical system in a chaotic region, where due 
to the non-linearity of the equation, very small initial perturbations are rapidly 
amplified and the distance between similar trajectories diverges exponentially. 

In our case this behavior occurs even if do = = 10~ 4 , see fig.4. In general 
one finds an increase of the populated wave numbers together with a strong 
mixing as the evolution proceeds, despite the fact that only one wave number k 
is initially occupied. The final frequencies spectrum has little resemblance with 
the initial one. Note however that the initial mode is not completely damped, 
and some memory of it is still retained. We have checked that this general trend 
depends neither on the initial nor on the shape of the initial profile. 

In fig. 3a it is displayed the time evolution of the distance between trajectories, 
as defined in eq.(4), normalized with respect to its initial value do = 10~ 4 . We 
note that, while the distance remains constant in the regular zone as time goes 
on, it increases exponentially when the system starts in the spinodal region. In 
fig. 3b \(t) as defined in eq.(2) is plotted vs. time. As one can see, A tends to 
zero when dynamics is regular, indicating that the trajectories actually do not 
diverge between each other. On the contrary, A rapidly converges to a value A 
which is not zero in the unstable regions and unambiguously defines the largest 
Lyapunov exponent A. We note that our time scale for the onset of chaos is in 
agreement with other numerical simulations on instability growth ]7], ^|, and 
with recent experimental results on multifragmentation [jnj. 

We have checked that the results are stable with respect to the increase of 
the numerical precision. This is shown in fig.4, where the time evolution of A is 
plotted for different values of the time integration step and initial distance. No 
appreciable variation of the results is obtained using a value of do smaller than 
10~ 4 or decreasing the time step. Therefore the divergence of the trajectories is 
an intrinsic dynamical property of the nuclear mean field in this region. 

The inclusion of two-body collisions, performed by solving the VN equation as 
done in ref. M, but without the Langevin term, does not modify substantially the 
dynamics, as also found independently in ref.0. This can be seen in fig. 5, where 
A is displayed vs. time for a density well inside the spinodal region. We note only 
a slight variation - about 10 %- of A when including the collision term. In fact, 
since the dynamics is chaotic, the exponential divergence between trajectories 
can hardly be affected by the collisions, which have a small effect at low density 
because of the diluteness of the gas ^|, |7j . 

Finally in fig. 6 we show the extracted values of A calculated at different den- 
sities. The values of A increase from zero at p = 0.66p to a a flat maximum 
around the density p = 0.35p - The curve reported in fig.6 , which is essentially 
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independent from the initial k, clearly defines the nuclear matter spinodal region 
up to its upper limit. 

It is worthwhile to compare the characteristic times tmf — h/\ which de- 
fines the time scale of the divergence between mean field trajectories, with the 
single particle characteristic time r sp = ti/Ep, being Ep the Fermi energy at the 
given density. We have found that, for densities p < 0Ap , the divergence time 
is smaller than the single particle time, and therefore in this region the notion 
itself of mean field ceases to have any validity. 



4 Conclusions 

In the framework of a purely mean field theory, we have shown that the nuclear 
dynamics inside the EOS spinodal region is chaotic. Without adding any noise, 
and even without two-body collisions, spontaneous symmetry breaking occurs 
leading to large unpredictable density fluctuations. We characterize the chaoticity 
of the system through the largest Lyapunov exponent, for which we suggest a 
computational recipe. The latter does not depend on the particular choice of the 
initial conditions or the time integration step, but it is strongly related to the 
initial density. The small influence due to two-body collisions has been discussed, 
as well as the validity limits of the mean field theory. These results confirm 
quantitatively a scenario already suggested by other authors || |, [7|, |TTJ and can 
have important consequences for the whole picture of nuclear multifragmentation. 

Generally speaking, a central collision between heavy ions at intermediate en- 
ergy can be viewed as follows. In the first stage of the reaction nuclear matter gets 
compressed and a highly excited composite nuclear system is formed. During this 
phase, the dynamics is insensitive to small modifications of the initial conditions. 
Once the maximum compression is reached, nuclear matter starts to expand and 
can merge into the spinodal region. According to the results presented in this 
paper, at this point the dynamics becomes chaotic and the corresponding large 
density fluctuations dominate the cluster formation. The degree of chaoticity and 
fluctuations depends of course on the details of the dynamics which brings the 
reaction inside the spinodal region. We would like to stress that, since initially 
very close events follow completely different dynamical paths, an event-by-event 
description is strongly advised. Averaging over different theoretical paths can 



partially mask the non-linear evolution [12]. This scenario can support the be- 



havior suggested by percolation model and at the same time disagree with an 
equilibrium statistical description since some memory of the initial mode seems 
to remain. Unfortunately, at the moment, the extension of these calculations to 
the case of collisions between heavy ions appears problematic even in two dimen- 
sions. A realistic theoretical simulation which is able to describe how the system 
enters into the spinodal region and how long it stays there is desirable. On the 
other hand an unambiguous experimental signature is needed. 
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Figure captions 



Fig.l Evolution of the density profile (a) and its Fourier transform (b) at initial 
normal density, see text. 

Fig. 2 Evolution of the density profile (a) and its Fourier transform (b) at initial 
density p = 0.5po, see text. 

Fig.2bis Evolution of the density profile (a) and its Fourier transform (b) at 
initial density p = 0.51po 5 see text. 

Fig. 3 For three values of initial density the evolution of the distance bet- 
ween two trajectories (a) and the corresponding X(t) (b) are shown. The initial 
distance is do = 1CT 4 , see text. 

Fig. 4 The quantity X(t) as a function of time changing do and the time step. 
Using a value of do smaller than 10~ 4 no change occurs. 

Fig. 5 Behavior of X(t) as a function of time with and without collision term 
in eq.(l). 

Fig. 6 The Lyapunov exponent A vs. p/p at temperature T=3 MeV, see text. 
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